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I. INTRODUCTION 


Accurate path control of ships and underwater vehicles along prescribed 
geographical paths is a fundamental problem which is becoming increasingly important, 
particularly as the missions of ocean vehicles become more sophisticated with strict 
requirements for performance. In order for a control law to be able to perform its 
mission in a realistic operational scenario it has to be robust enough so that it can 
maintain stability and accuracy of operations in the presence of modeling errors and 
environmental uncertainties. The robustness properties of the design are particularly 
important due to the unpredictable nature of the ocean environment and the changes in 
the hydrodynamic characteristics of the vehicle during turning, changes in the forward 
speed, or operations in proximity to other objects in the area. For these reasons, there 
exists a need for the analysis of the robustness characteristics of a particular control law 
design and the establishment of a rational operational envelope based on stability and 
performance criteria. Previous studies (Parsons & Cuong, 1977) showed that gain 
adaptation is highly desirable due to changes in the linearized vehicle hydrodynamics 
with different operating conditions, such as depth under keel. The resulting adaptation 
scheme (Parsons & Cuong, 1980) required significant vehicle motion which may be 
undesirable when operating in restricted water, or in object recognition and localization 
tasks. Integral control techniques (Parsons & Cuong, 1981) proved quite effective, but 


neglected the behavior of the vehicle which becomes very important at low speeds and 


hover operations. Model based compensators exhibit robust behavior under conditions 
of parameter uncertainty which is as good as the classical linear quadratic regulators for 
linear output feedback systems (Healey, 1992). Alternatively, sliding mode controllers 
exhibit very robust characteristics given an estimate of the parameter uncertainty and/or 
disturbances (Papoulias & Healey, 1992), (Yoeger & Slotine, 1985). Sliding mode 
control, however, does not offer an infinitely robust design, and it suffers from a series 
of bifurcation phenomena and loss of stability unless proper care is exercised (Papoulias, 
1991). 

In this work we analyze the problem of loss of stability of a path keeping control 
law under conditions far from nominal. We assume that the autopilot has been designed 
based on a nominal model, whereas the actual system is different. For demonstration 
purposes we employ a linear full state feedback control law, but the methods are quite 
general and can be used for other designs as well. The main loss of stability cases 
analyzed We occur in the form of generic bifurcations to periodic solutions 
(Guckenheimer & Holmes, 1983). We use center manifold reduction techniques and 
averaging in order to capture the stability properties of the resulting limit cycles (Chow 
& Mallet-Paret, 1977). Particular emphasis is placed on the control gains as the primary 
bifurcation parameters, since they are related to gain margins in linear control theory 
(Friedland, 1986). We aT extensive use of numerical integrations in order to confirm 
the theoretical results. All computations in this work are conducted for the NPS 


autonomous underwater vehicle (Bahrke, 1992), and all results are presented in standard 


dimensionless quantities with respect to vehicle length, 7.3 ft, and nominal forward 


speed, 2 ft/sec. 


Hi. PROBLEM FORMULATION 


A. INTRODUCTION 

The equations of motion of an ocean vehicle in the horizontal plane are presented 
in this chapter. A linear feedback control law is designed based on the linearized 
equations in yaw, sway, and rate of change of heading angle and lateral deviation error 
to provide path keeping. Loss of stability is examined for small changes from nominal 


in feedback gains, and system properties. 


B. EQUATIONS OF MOTION 
The mathematical modeling of a steering system of a vehicle consists of the 
nonlinear sway and yaw equations of motion. Newton’s equations of motion in a moving 


coordinate frame fixed at the ships geometrical center are’ 


m(v+r+x fF) =YF+Y vt r+YvrV¥, 6 é +Y, 5 P 
Stern 


- f CoACCv+ Er) |v+Cr| ae 


bow 


[r+mx(v+r)=N,+Nv+N +N v+N, 6.4+N 3,0 ; 
bow 9] 


- [ Cyh(C)(v+ Cr) |v+Cr|dv 


Stern 


' Equations are nondimensionalized with constant forward speed u,ship length L, and the dimensionless time 
being t.u/L. 


where only nonzero terms are kept in the model and the symbols are explained in the 
nomenclature. The cross flow integral drag terms in the equations of motion are very 
important in low speeds, but for higher speeds their effect is much smaller and the 
steering response predominantly linear. Also for maximum maneuverability, the vehicle 
bow rudder is deflected at the same amount and opposite to the stern rudder. 


Equations (1), and (2) can be written in the form 


V=a, ,V+a,,r+b,5 (3) 


F=a,,V+a,,r+b,6 (4) 


where 


Da,,= (,-N,)¥Y,-(mx,-Y,)N, 

Da,,= (1,-N,)(Qn-Y,)-(mxg-y,)(n,-mXq)_ , 

Da,,= (m-Y,)N,-(mx,-N,)Y, 

Da,,= (m~Y,)(N,-mx,)-(mx,-N,)(Y,-m), (5) 
Db, = Cy JY 3, 3) —(mx,-Y,)(N; -N, ) ; 

Db, = (m “YN, -N 3) > 

D (,-N,)(m-Y,)-(mx,-Y,)(N,-mXq)_ . 


The transfer function between rudder angle 6 and yaw angular velocity r is obtained from 


Equations (3) and (4). 


£_bSMayby-ayb) 6 
a) S2 -(a,, +a,,)S +(a, 1222 ~@4,1) 


For low frequency maneuvering motions this second order equation can be 


approximated by expanding in Taylor series and keeping only the first order terms only 


or fF = ar+bé (7) 


where 


(2, ,4,.-@,.4,,)(a+21b1 -a,,b,) 
(4, +45y)(4,,0, -€, ,b,) +b, (4, 45. - 4441) (8) 
(a,,5, Ele 


(, +255)(@,,5, -@, 5, -4, b,) +b, (A, ;25y -A,545,) 


Equation (7), Nomoto’s Equation is very useful in control system design since no 
sway velocity feedback necessary, and it gives the fundamental turing performance. At 


steady state, Equation (6) turns out to be 


r=0Q 
ar+b6 = 0 (9) 
r = a5 

- Qa 


From Equation (9) it can be seen that the relationship between yaw rate and rudder 
control angle are linear with a slope -b/a. Experimental results show that the actual 
relationship is not linear. Large increases in rudder angle fail to increase the yaw rate 
according to the amount predicted from Equation (9). 

The linear change in the time derivative of yaw rate due to rudder angle can be 
augmented by a nonlinear term of the form ’a,r°’. This is introduced to model the 
appropriate speed loss during turning. The term ’a,’ has the same sign as ’a’. The new 


yaw equation is therefore 


Fear+a,r>+bd (10) 
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Figure 1 Difference between Nomoto’s equations and experimental results. 


The model becomes complete with the expression of yaw rate and the inertial 
deviation rate from the commanded path. Therefore the coupled nonlinear equations of 


motion for the marine vehicle in the horizontal plane are 


y = (11) 
fF = ar+a,r>+bé (12) 


This system of equations forms the basis for the controller design. 





Figure 2 


C. CONTROL LAW 
Equations (11), (12), (13) govern the steering control of the model used in this 
section.The control law can be expressed as, 


6 = 8 tan") (14) 


Sat 


where in the vicinity of » = r = y = O we have 


5, = y¥ +K r+K,y (15) 
6 is the rudder angle and Ky ,K,, K, are the control gains of the system. The linear 
control law is 6, . The rudder angle 6 is defined by a hyperbolic tangent function to 
include the saturation to our problem as shown in Figure (3). Saturation occurs at 6,,, 
which is the saturation limit generally can be taken as 0.4 rad. The linearized form of 


equation of motions in the vicinity of y = r = y= 0 are 


ber (16) 
r =ar+bd, (17) 
yay ee 





Figure 3 


These equations can be expressed in space state form as 


X = AX (19) 


A is the Jacobian matrix of the system and X is the state vector 


0 1 


0 
A = |bK, a+bK, bK, (20) 
1 0 0 
X=[Wrny]? = 
The characteristic equation of the matrix A is 
A? -(a+bK,)A’ -bK, 2.-bK, =0 (22) 
The desired characteristic equation has the form 


AP +07 +0, A +045=0 (23) 
The gains can be generated in terms of the coefficient of actual and desired characteristic 
equation from Equations (22) and (23), 





ae (24) 
Vv b 

K - _ ta (25) 

j b 

ge ge 


y 


(26) 


The desired characteristic equation can be written with respect to the desired 


natural frequency and some optimum coefficients. The ITAE criterion for a third order 
equation is (Dorf, 1992) 
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$?+1.75_ S?+2.150,S +, =0 (27) 
Therefore the control gains can be calculated for a given natural frequency, as 


2 


a) = @, 
w= 2.1502 (28) 
a, = 1.750, 


D. LOSS OF STABILITY 


1. Introduction 
The control law guarantees stability if all real parts of eigenvalues of the 
jacobian matrix are negative. Equation (23) is stable with the chosen control gains. A 
small perturbation in the control gains or coefficients of the system affects its stability. 
In this section we presented the computations for the critical value for change in the 


gains and the coefficients of the system. 


2. Stability Considerations 
For a system with a third order characteristic equation, the Routh Hurwitz 


criterion (Dorf, 1992) requires for stability 
@, &, =O, (29) 


a. Perturbation in K, 
The change in K, is_ defined by a coefficient c, where c is any real 


number. The linearized equations of the system for a small change in K, are 


11 


W =r 
* =ar+bd, (30) 
ae 
where the control law is 
6, = cK, w+K r+K yy (31) 
Therefore the characteristic equation is 
A? -(a +bK,) M*-cbK,A -bK, = 0 (32) 


Equation (32) should meet the Routh Hurwitz criterion for stability. 


(a+bK,)cbK, 2 -bK, (33) 


The critical value of c where a_ stability change occurs in the system is presented in 


Equation (34) 





Cae = et (34) 
TH, Xa, 
Using the coefficients, the critical c is 
C crite, = 0.2658 (35) 
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Gro52 


stable region : 


0.266 


Ccritical 


unstoble region 





wn 


Figure 4 C.,, WS natural frequency for Ky. 


b. Perturbation In K, 


Similarly for perturbation in K, the characteristic equation is 


A? -(a+cbK,)d’ -bK, A -bK, = 0 (36) 
The critical value of c, 
ee Ol (37) 
E, — (5 +a) 
or 
@,+2.135a (38) 
ees 
eee Ol fo. 4a) 
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stoble region 


Ccritical 





Oo 0.5 1 ie 2 2.0 3 


wr 


Figure 5 c,,;, versus natural frequency for K,. 


c. Perturbation In K, 


The results for changes in K, are presented similarly. The characteristic 
equation, is 


\>-(a+bK,)?-bK, 4 -bcK, = 0 





(39) 
The critical value of c is 
_ &) &, 40 
C. rt, = " (40) 
Ccrity = 3.7625 (41) 
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525 


unstable: region 


o2765 


Ceritical 


stable région 
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wn 


Figure 6 c,,, versus natural frequency for K,. 


d. Perturbation in a 


After similar calculations the critical value of c is 


(0, +2) %,-Gp 


C rit a a (42) 
2 ad, 
Or 
28484 
Soe Tae (43) 
e Qa 


13: 


Ccritical 





wn 


Figure 7 ¢,,,, versus natural frequency for a 


e. Perturbation in b 


After similar calculations as with the other cases the characteristic 


equation 1s, 


NP +(c(a,+a)-a)a*+ca,A+ca,=0 (44) 


and the critical values are 


ie (45) 
a, (a,+a) 


C= 
cril, 


Or 
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w +2.15a 
C crit = SSS (46) 
eee 19( 1-75 © +a) 
Examining Equation (46), we note that all the coefficients of polynomial 


must have the same sign, if all the roots are in the left hand plane.The necessary 


conditions for stability are, 


c>0 
(47) 


c(a,+a)-az0 


The system is unstable for all values of the natural frequency if c is less then zero. When 
the natural frequency is in the interval 0 < w, < |a/1.75] stability changes from stable 
to unstable with an increase in c (Figure (8)) for positive values of c. In the interval 
w, 2 |2.15a| stability changes from unstable to stable while c crosses the critical point 
again for positive values of c. In both cases a pair of complex conjugate eigenvalues 
crosses the imaginary axis. Stability is also lost when c becomes negative, but this is 
associated with a real eigenvalue crossing zero. In this work we concentrate on the c>0 


Case since in applications it is unlikely that a change in the sign of b will occur. 


ie 


Ccritical 





fe) 0.5 1 Wee 


2 ANS - 3-9 


Figure 8 c,,, versus natural frequency for b. The c,.., goes to infinity asymptotically at 
w, =-a/1.75. 
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I. HOPF BIFURCATION 


A. INTRODUCTION 

Hopf bifurcation is the simplest bifurcation, in which under the variation of a single 
control parameter a stable focus equilibrium bifurcates into an unstable focus surrounded 
by a growing limit cycle. As c crosses the critical value one pair of complex conjugate 
eigenvalues of the linear system matrix crosses transversely the imaginary axis. In this 
generic Poincare’ - Andronov - Hopf bifurcation a family of periodic solutions coexists 
with the stable / unstable nominal equilibrium state. Locally as ¢ approaches c,,, the 
periodic solutions are located on the two dimensional Euclidean plane spanned by the 
eigenvectors of Jacobian matrix of the system which corresponds to the critical pair of 
eigenvalues. In this chapter stability 


properties of the periodic solutions are 


established. In order to establish those | r y 
properties the main nonlinear terms that 

dominate the system are isolated. Center O YY, 
manifold theory is used to reduce the flow 

to a two dimensional manifold. The Figure 9 


method of averaging is applied to the 
reduced system. Hopf bifurcations are examined for small changes in each gain and 


system dynamics coefficients. 
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B. THIRD ORDER EXPANSION OF THE SYSTEM EQUATIONS 


1. Perturbation in Ky 


In this case the equations of motions are 


Ww =r 
f= arta,r>+bd (48) 
y = siny 


where 


6 
5 = 5 sat : 
8 at (49) 


by = CK, +K r+ Ky 


The trivial equilibrium state is characterized by y=r=y=O .Taylor expansion of the 


nonlinear terms about the equilibrium gives 








ye 
siny = YF +O(5) 
5 3 (50) 
8 = 8,-——-+0(5) 
36.0, 
The system equations are 
Wer 
3 
; sap, - Bee 
= ar+a,r = 
bee all 
yw 6 
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In state space form they are written as 


X = AX+g7(X) (52) 


where A is the linearized system matrix, and g°(x) contains the nonlinear terms, 





0 1 0 
A=\cbK, a+bK, bK, (53) 
1 ac 
gi(x) = 0 
ame) 10M 
82%) = ~~ 5% (54) 
sat 
3 
g(x) = - 


where the 6,° term is computed as 


89 = cKy W°+K,7r>+K,y° 
307K, (K,r+K,y)+3K, 1? (cK, p +K,y) (55) 
3K, y*(cK, +K,r) +6cK, KK, ry 


With the computed gains from equations (24), (25), (26), Jacobian matrix is 


0 ] 0 
] 0 O 


The eigenvalues of matrix A are computed at the bifurcation point c,,, from equation 


(34) where a pair of complex conjugate roots with real parts are obtained. 
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7 a a 
1 
a 

2 

. (57) 
A, = -j |— 

a, 
a = -@, 


For the above system of equations a transformation matrmx of eigenvectors can be 


introduced, 





1 0 1 
o 
ies - 5 (58) 
0 |— -Ia, 
| Xo 
The transformation 
x=T2 z=T 'x (59) 


transforms the system into its normal coordinate form, 


Z=T Alz tie) (60) 


where, 


ZZ 





a. \7/2 
jiel =) 2) ere 0 
a (a3, 
0 
and 
jae |e 
Oo 
A= |i 
ae 0 
Oo 
0 0 


The transformation of physical variables to normal coordinates is 


x= Tz, 


Y = 2, +2, 


a, 
a 


us 





” (61) 


(62) 


(63) 


Y 0 Tite A <1 
ig = a, : 4 (64) 
y 


; Hs PML Ss 
Xo o> 


The coordinate z, corresponds to the eigenvalue A; which is real and negative. 





A,, A, denote the complex conjugate pair with zero real parts. At the critical point two 
distinct eigenspaces are spanned by the eigenvectors associated with the two sets of 
eigenvalues. These eigenspaces can be viewed as local approximations of the invariant 
manifolds the center and the stable manifold. The system exhibits its essential 
bifurcational behavior on the centre manifold associated with A, and )\,. Centre manifold 
theory reduces the flow to a two dimensional manifold (Hopf bifurcation). According to 
centre manifold theory the coordinate z,; is expressed in terms of z, and z,, but this 


expression is higher order. Equation (64) can be written as, 


Waa, 
6 Ng, 2 (65) 
» ae 

aa 


Substituting the equations (65) to (60) the system equations in normal coordinates 


become 
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a 


ae oO, | a es 2 3 
4; = FA Zp 92 St 432 2n tM gk 
z (66) 
; Xo : an Do ae 3 
= a Lyt 4321 119921 & 119341 <2 Tyg ho 
2 


The terms r,; are computed later from our equations at the bifurcation point c,,;,. For 


values of c in the neighborhood of the critical point, equation (66) becomes, 


7 = B/ =| atte “o + ; + Z + ae ‘ 
Z, = Bez,-| we SH Sy FP 9h tt 132122 1 4k 
(67) 
a vere “0 rae! + un cae ar : 
~% = |e a Z,+ P°EZ, +741 Zy +0 y9Zy Z +1932 2 +0 y42 


74 


The parameter e is the difference between c,,,, and c or, c=c,,,, +¢. The terms $' and w' 


crit 
are the derivative of the real and imaginary parts of eigenvalues with respect to c 


evaluated at c=c,,,,, and computed from a perturbation series approach. The perturbation 


is expressed as follows from equation (34), 


Xo 


c= +€ (68) 





ae Ve 
The characteristic equation can be written as follows for a small change for K, in the 


neighborhood of c...,. 


o 
M3 +0,A7+(—+0,€)A+a, = 0 (69) 
o 
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The solutions of equation (69) are presented in equation (70) 


Aig = beF/J “Oe 
0. (70) 


A, = -@,+a,€ 
The variables b, and c, can be found by substituting equation (70) to (69) and neglecting 
the successive powers e’,e°. 


2 


| 
ys oe 
2(a,+%>) 
(71) 

ed 

4%, | 

de} 

a a a Oe 
2( +>) 


6' and w' can be expressed by the definition of the derivative where « is zero in the limit, 





b,e-0 
Cue = 
€ 
OL OL (72) 
—+¢,e|- _o 
OL OL 
wo! = tL tf Nt xg 


Equation (72) can be written in terms of the desired natural frequency, 


TO 
— 
U 


-0.6841 w, 
|= SOG 


(73) 


€ 
T 
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a. Calculations of r,, Terms 


The terms 1; are generated by the last term of equation (60), where 


T'g*(Tz) is the third order nonlinear part. From equations (61) and (54) 


eae ee ee) 2s | 0 0 
0.5 (a3 +0) Xo ae) ; 3/2 ? 
3/2/,,3 0 3 3 
a (a+ 1 -aJ—| z-——6 
rem Aes] 9 tg, || AB] Ee a 
a5 2 
0 ier jac Be 
Xo Xs : 


where, 


Pa psi ee ae 
i, = ric +a 2] orie,-3f 3 ote rail 2) z} (75) 
b? a. 4 a 


2 2 a. 


After expanding the first row of equation (74) and organizing the terms, we get 





OL a, \* 
r..= —2 (a5 + 0,) SU ae tt (76) 
ea Spa N%,) 6 
sat 
OE 
SO. & {3 (05 + cy) = (77) 
v4 252 
2 b* 8. 


a’ 


3 
Ty, = (@, +0) 





(78) 


§ 722! 


Zi 





0 
bie {Z| arenes (79) 


0 35287 


Sat 


Similarly after expanding the second row of equation (74) and organizing, we get 











; 3? Tone 
p= -(0. 0) ———— rr eee (80) 
21 200 Se s2. seu 
3b7 5. %5 ho 
o 
2s (09 + a9) — (81) 
Z b Sat 
1/2 
ao 3 ae 
To, = Ser) a (82) 
06.5 b*6_., 
05 +0 a? 
1, = NG (83) 
24 3 9 2 
C5 3b°8_, 
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b. Averaging 


Equation (67) can be written in a simpler form , 


ot 
Z, = B’ez,-| w’e+ E 
cL 
2 
; 'e+ |—o| 24 Blez, +F ) 
Oye e RS ae =a €Z, + F, (2,2, 
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where F, and F, contain the third order expansion terms. The use of polar coordinates 


2, + F,(Z,,2)) 





(34) 


makes it possible to decouple equation (84) . We use the polar transformation 
Z,= R cos(@) . Z,=R sin(@) 


After some algebra, equation (83) can be expressed in terms of R and @, 


R = B’eR+F,(Rcos(6), Rsin(6))cos(6) + F, (Reos(8), Rsin(6))sin(@) (85) 


(36) 
R6 = {it-0 R +F,,(Rcos(6), Rsin(6))cos(6) - F (Rcos(6), Rsin(6)) sin(8) 


Equation (85) can be written as follows, 


R = B/eR+P(0)R? (37) 


where P(@) is 2a periodic in the angular coordinate @, 


P(®) = r,,cos*(8) +r,,cos(®)sin(0) +r,,cos?(6)sin?(8) +r, ,cos(®)sin3(6) 
+r, ,cos’(8)sin(8) +7,,cos"(8)sin7(6) +7,,cos(6)sin*(6) +7,,sin*(8) 


Equation (87) is averaged over one cycle to obtain an equation with constant coefficients. 


R = B/eR+KR? 


where K is defined by, 
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Equation (90) is simplified as equation (91) after evaluation of the integral. 
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Equation (91) can be expressed in terms of the natural frequency, 
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Equation (87) can be expressed as, 
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where, 


F(8) = r,,cos*(6) +7r,,cos*(6)sin(@) +r,,cos*(8)sin*(6) +r,,cos(6)sin*(8) 
-r, ,cos*(6) sin(6) -7,, cos?(8) sin?(6) - r,,cos(6) sin*(6) - r,,sin*(6) 


After averaging equation (93), a constant coefficient M 1s formed, 
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where M is defined by, 
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2. Perturbation in K, 
The calculation of the formulas for the parameters K,, K,, a and b is similar 
to the calculation for K,. Therefore in the following sections only the main results are 


presented. The system equations for a change in K, are 
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y = sinpy 
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The Jacobian matrix is 


ak 


0 1 0 
A =\-e, (Cae conl a, 


1 0 0 
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The eigenvalues of the Jacobian matrix are calculated at the bifurcation point from 


equation (38), 


t (101) 


The state vector is presented in normal coordinates where z, corresponds to Az, 
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r= -/a, % (102) 
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The system equation in normal coordinates are, 
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The coefficients of the third order expansion terms are, 
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The cubic coefficients K and M are obtained after averaging , 
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3. Perturbation in K, 


The formulas for perturbations in K, are presented in this section. The 


equations of motion are, 


Wer 
fF = ar+a,r°+bd (115) 
y = sin 
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The Jacobian matrix, 
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The normalized equations of motion, 
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The coefficients for the third order terms, 
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The cubic coefficients K and M , 
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4. Perturbation in a 


The equations of motion, 


w=r 
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y = sinp 


where, 
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The Jacobian matrix, 
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The eigenvalues at the bifurcation point, 


Ais = FJ /e, 
OL 

——— 
Oy 


and the state equations in normalized coordinates, 
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The normalized equations of motion in the neighborhood of the bifurcation point, 
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The cubic coefficients K and M are 
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5. Perturbation in b 


The equations of motion are, 
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The cubic coefficients K and M are expressed in terms of the natural frequency w, 
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C. RESULTS 

Before we present the results of analytic study of Hopf bifurcation we have to 
assess the effect of the a, term. The a, term plays a role in equation (12), and from the 
nature of the equation we can see that a, has to have a negative value. This is because 
of the softening spring characteristic of the steady state r - 6 curve, as shown in 


Figure 10. 





Figure 10 r -6 graph for different a, values. 


In the previous sections expressions for 6’ and K were found for equation (89). 
The value of K depends on only the nonlinear terms and since the eigenvalues of the 
Jacobian matrix cross the imaginary axis with nonzero speed the term £° is nonzero. 
Equation (89) has two steady state solutions, one at R=0O which corresponds to the 


trivial equilibrium solution at zero and one at 


r2--£¢ (168) 


This equilibrium solution corresponds to a periodic solution (limit cycle) in the cartesian 
coordinates z,, Z,. From equation (168) we can conclude that, 
1. If 8° >0 then, 
a. If K>0, then unstable period solutions coexist with the stable equilibrium 
for «<Q, and 
b. If K <Q, then stable period solutions coexist with the unstable equilibnum 
for e>0, 
2. If 8° <O then, 
a. If K>0, then unstable period solutions coexist with the stable equilibrium 
for e>Q, and 
b. If K <0, then stable period solutions coexist with the unstable equilibrium 
for «<Q. 
The stable periodic solutions form the supercritical Poincare” - Andronov - Hopf 
(PAH) bifurcation, while the unstable periodic solutions form the subcritical PAH 
bifurcation. The period of the limit cycles is computed by substituting equation (168) in 


(94) 
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Figure 11 The two general types of Hopf bifurcation. 
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The existence and stability of periodic solutions in our cases is examined in the 
following subsections. The graphs of w° and @° are presented in figures (12) and (13). 
We note that for b and for natural frequency in the range |a/1.75| <w,<|2.15a| the 
system does not have a complex pair of eigenvalues. Hopf bifurcation does not occur in 


that interval. 
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Figure 12 6° versus w,. 





Figure 13 w° versus ),. 
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1. Perturbation in K, 

The graph of 6° for Ky is shown in figure (13). 6° is always negative for K,, 
and the solutions of equation (91) is presented in figures (14), (15), (16). K is less than 
zero for low natural frequencies. For those frequencies we get a supercritical Hopf 
bifurcation and a stable periodic solution exists. However when the natural frequency 
increases the bifurcation shifts to a subcritical Hopf bifurcation. Figure (14) is obtained 
for a; = 0 and various 6,,,. Figure (15) is obtained for a, = -3 and vanes Ost Figure 
(16) is obtained for 6,,,=0.4 and various a,;.When |a,| increases the domain of the 


natural frequencies, which causes supercritical Hopf bifurcation, increases as well. 





Wn 


Figure 14 K,, versus w, for a,;= 0 and various 6,,. 
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Figure 16 K,, versus w, for 6,,,=0.4 and various a3. 
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2. Perturbation in K, 

The term 6° for K, changes its sign at w,= |a/1.75| and can be observed 
from Figure (12). Therefore, for values of K<0O _ there exists a supercritical Hopf 
bifurcation whereas subcritical Hopf bifurcation forms when K changes its sign. From 
Figure (17) it is observed that 6,,, does not affect the domain of the natural frequency 
for supercritical Hopf bifurcation. However, an increase in the |a,| term increases the 


domain of supercritical Hopf bifurcation, as can be seen from Figures (18) and (19). 





Figure 17 K,, versus w, for a,;=O and various 6,,. 
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Figure 18 K,, versus w, for a,=-3 and various 6,,. 
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Figure 19 K,, versus w, for 6,,,=0.4 and various a;. 
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3. Perturbation in K, 

For this case 6° is always positive and for ¢€>0O there exists an unstable 
equilibrium, therefore for K<O a supercritical Hopf bifurcation exists. Similarly for 
K >0a subcritical Hopf bifurcation exists. An increase in 6,,, increases the domain of the 
supercritical Hopf bifurcation over the natural frequency, see Figure (20). Also an 
increase in |a,| increases the domain of the supercritical Hopf bifurcation, see Figures 


(21), (22). 
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Figure 20 Ky, versus w, for a,;= 0 and various 6,,, 
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Figure 21 K,, versus w, for a,= -3 and various 6,,, 





Figure 22 Kx, versus w, for 6,,,=0.4 and various a3. 


53 


4. Perturbation in a 
6° is always positive for this case. For «>0 and K <0 there exists a stable 
period solution coexisting with the unstable equilibrium. When a,= 0, a change in 6,,, 
does not affect the supercritical Hopf bifurcation region of the natural frequency, Figure 
(23). If a; is non zero, an increase in 6,,, increases the domain of the natural frequency 
which causes the supercritical Hopf bifurcation. An increase in |a,| also increases the 
domain of the natural frequencies for supercritical Hopf bifurcation, see Figures (24), 


(25). 





Wn 


Figure 23K, versus w, for a,;= 0 and various 6,, 
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Figure 24 = K, versus w, for a;= -3 and various 6,,, 





Figure 25K, versus w, for 6,,,.=0.4 and various a. 
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5. Perturbation in b 

Previously it was shown that there is no Hopf bifurcation in the interval, 
|1.75/a| <w,<|2.15a] K can be examined in two regions, 0<w,<|1.75/a| and 
w, > |21.15a|. In the first region 6“ is always negative, «>0, and K is always negative. 
In this region only supercritical Hopf bifurcation occurs. Changes in 6,,, and |a,| have 
no effect in this region, see Figures (26), (27), (28) . For the second region B “ <0, «<Q, 
K <0 initially and supercritical Hopf bifurcations exist. For higher natural frequencies, 
subcritical Hopf bifurcation is observed. Figures (29) and (30) show that changes in 6,,, 


or a; have no significant effect. 
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=().4 and various a; for the first region. 
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Figure 26 K, versus w, for 6 
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Figure 27 K, versus w, for a,= 0 and various 6,,, for the first region. 
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Figure 28 K, versus w, for a,= -3 and various 6,,, for the first region. 
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Figure 29 K, versus w, for 6,, =0.4 and various a; for the second region. 
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Figure 30 K, versus w, for 6,,= 0.4 and various a; for the second region. 
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IV. SIMULATIONS 


To illustrate the dynamics of the Hopf bifurcation we simulated the case in which 
there is a change in the gain K,. We compared the simulations with the analytic studies 
presented in previous chapters. In our case, loss of stability occurs at c=3.7625 from 
stable to unstable while c increases its value, as can be observed from Figure (6). In the 
following simulations we use the typical saturation angle for the rudder of 0.4 radians. 
The a, term is assumed to be zero. 

K is less than zero for the natural frequency 0.4 in Figure (20), therefore a 
supercritical Hopf bifurcation exists in this region. The behavior of this supercritical 
Hopf bifurcation is simulated in Figures (31), and (32). In Figure (31) cx,=3.5 is 
simulated for two different initial lateral deviations ( y) = 0.05 and y, = 0.5). 
Regardless of the initial conditions both curves converge to zero for this value of c. 
There is a unique stable steady state solution for this value of c. Figure (32) is simulated 
for K < 0 and w,= 0.4 and the same initial condition of y>=0.5. When c= 3.5 the 
simulation converges to zero as explained in Figure (31). When c is increased above the 
critical value (here c=4.0) we observe that the solution converges to a limit cycle. This 
stable limit cycle coexists with the unstable equilibrium solution. 

Figures (33) and (34) are presented to show the effects of a subcritical Hopf 
bifurcation. For this purpose w,=1.2 is simulated for Figures (33), and (34). w,=1.2 


is in the region where K>0O in Figure (20). The initial condition for this case is 


Se, 


Yo= 0.05. When c is set to 3.5 the solution converges to zero, but for the case where 
cis set to 4.0 the simulations converge to a stable limit cycle. This is the case in Figure 
(11-b). The simulation misses the unstable limit cycles and converges to the outer limit 
cycle shown in Figure (11-b). In Figure (34) the effects of the initial conditions are 
further illustrated. In Figure (34) for y,>=0.05, the simulations converge to zero, but for 
Y)=0.5 the simulations converge to a limit cycle. For this simulation the initial condition 
was large enough to shift the results to a limit cycle rather than the zero steady state 
solution. 

From the simulations we observe that in subcritical bifurcations the magnitude of 
the stable limit cycle is greater than the magnitude of the limit cycle in supercritical 
bifurcations. Similar results can be obtained for the other cases of variation in parameters 


a, b, and gains K, and K, . 
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SUPERCRITICAL HOPF BIFURCATION 
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Figure 31 Supercritical Hopf bifurcation for w,=0.4 , c=3.5 for two simulations 
with yg=0.05 and yo=0.5. 
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Figure 32 Supercritical Hopf bifurcation for w,=0.4 , yo= 0.5 for two simulations 
c=3.5 (stable equilibrium ) and c=4 (unstable equilibrium) . 
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SUBCRITICAL HOPF BIFURCATION 





=o Se ee ae 





0 10 20 30 40 350 60 70 80 g0 100 
time 


Figure 33 Subcritical Hopf bifurcation for w,=1.2 , yp= 0.5 for two simulations 
c=3.5 (stable equilibrium ) and c=4 (unstable equilibrium) . 
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Figure 34 Subcritical Hopf bifurcation for w,=1.2 (K >0) where c=3.5. The simulations 
for y)=0.05 converged to the steady state value and for y,)=0.5 to the limit cycle. 
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V. SUMMARY AND CONCLUSIONS 


An analytic investigation of the nonlinear dynamic response characteristics of a 
steering control law of marine vehicles has been presented. Bifurcation theory techniques 
were utilized in order to assess the behavior of the system upon initial loss of stability. 
The main bifurcation parameters were the natural frequency, control saturation level, 
system parameters and gains. The main conclusions of this research can be summarized 
as follows: 

1. There exists a critical point for each gain and system parameters for stability of 
straight line motion. The loss of stability occurs in the form of Hopf bifurcation for each 
case. As the parameter crosses its critical value, a family of periodic orbits ("self 
oscillations") develops. 

2. The desired characteristic equation from the ITAE criterion has a faster response 
for higher natural frequency. However, subcritical Hopf bifurcation develops when the 
natural frequency is sufficiently high. If ¢ is increased past c,,; the amplitude of lateral 
deviation undergoes a jump and takes "large" values ( "hard loss of stability or hard 
generation of limit cycles"). For lower natural frequencies, stability loss occurs as 
supercritical where the limit cycle amplitude grows continuously (" soft loss of 
stability”). 

3. The stationary point where change of stability occurs is related only to the natural 


frequency of the system and the hydrodynamic coefficient a. 
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4. Higher saturation limit of the rudder angle increases the range for supercritical 
Hopf bifurcation. Theoretically, if a saturation limit does not exist all bifurcations 
develop as supercritical. 
5. The term a; has a similar effect as the saturation limit of rudder angle. An increase 
in the |a,| term increases the natural frequency where supercritical Hopf bifurcation 
exists, but an increase in the |a,| term is not a desired condition since vehicle response 
to rudder angle slows down. 
6. Transitions from supercritical to subcritical are very important for the design of the 
control law. A supercritical Hopf bifurcation is preferred over a subcritical Hopf 
bifurcation since the subcritical Hopf bifurcation may develop a rapid dynamic jump to 
a new limit cycle. 

Some recommendations for further research are as follows: 
1. Comparative studies must be performed including the observer dynamics. The 
effect of reren in gains and system parameters must be studied to predict the 
accuracy of the observer. 
2. Continuation techniques for periodic solutions must be performed for further 
stability analysis. The cases for simultaneous variations of more than one terms must be 


studied to understand the dynamic response of the system including any nonlinear terms. 


64 


LIST OF REFERENCES 


Bahrke,K. On-line Identification of the Steering and Diving Response Parameters of an 
Autonomous Underwater Vehicle from Experimental Data, Mechanical Engineer’s Thesis, 
Naval Postgraduate School, Monterey, California, 1992. 


Chow, S. N. and Mallet-Paret, J. "Integral Averaging and Bifurcation," Journal of 
Differential Equations, Vol. 26, pp. 112-159, 1977. 


Crane, C. L., Eda, H., and Landsburg, A., Principles of Naval Architecture (ed. E. V. 
Lewis), The Society of Naval Architects and Marine Engineers, New York, 1989. 


Friedland, B., Control System Design, An Introduction to State Space Methods, Mc- 
Graw Hill, New York, 1986. 


Guckenheimer, J. and Holmes P., Nonlinear Oscillations, Dynamical Systems, and 
Bifurcations of Vector Fields. Applied Mathematical Sciences 42, Springer-Verlag, New 
York, 1983. 


Hassard, B. and Wan, Y.H., "Bifurcation Formulae Derived from Center Manifold 
Theory," Journal of Mathematical Analysis and Applications, Vol. 63, pp. 297-312, 
1978. 


Healey, A. J., "Model-based Maneuvering Controls for Autonomous Underwater 
Vehicles," Journal of Dynamic Systems, Measurement, and Control, Transactions of the 
ASME, Vol. 114, pp. 614-622, 1992. 


Papoulias, F. A. and Healey, A. J., "Path Control of Surface Ships Using Sliding 
Modes," Journal of Ship Research, Vol. 36, p. 2, 1992. 


Papoulias, F. A., "Bifurcation Analysis of Line of Sight Vehicle Guidance Using Sliding 
Modes,” International Journal of Bifurcation and Chaos, Vol. 1, p. 4, 1991. 


Parsons, M. G. and Cuong, H. T., Optimal Stochastic Path Control of Ships in Shallow 
Water, Office of Naval Research Report No. ONR-CR-215-249-2F, 1977. 


Parsons, M. G. and Cuong, H. T., Adaptive Path Control of Surface Ships in Restricted 


Waters, Department of Naval Architecture and Marine Engineering, Report No. 211, The 
University of Michigan, Ann Arbor, 1980. 


65 


Parsons, M. G. and Cuong, H. T., Surface Ship Path Control of Surface Ships in 
Restricted Waters, Department of Naval Architecture and Marine Engineering, Report 
No. 233, The University of Michigan, Ann Arbor, 1980. 


Yoerger, D. R. and Slotine, J. J., “Robust Trajectory Control of Underwater Vehicles," 
IEEE Journal of Oceanic Engineering, Voli. 10, p. 4, 1985. 


Seydel, R., From Equilibrium To Chaos, Elsevier, New York, 1986. 


Papoulias, F. A., "On the Nonlinear Dynamics of Pursuit Guidance For Marine 
Vehicles," submitted for publication to Journal of Ship Research, 1992. 


Carr, J., Applications of CENTER Manifold Theory. Applied Mathematical Sciences 35, 
Springer-Verlag, New York, 1981. 


Marsen, J. E. and Mc Cracken M., The Hopf Bifurcation and Its Applications,. Applied 
Mathematical Sciences 19, Springer-Verlag, New York, 1976. 


Thompson, J. M. T. and Stewart, H. B., Nonlinear Dynamics And Chaos, John Wiley 
and Sons, New York, 1986. 


Dorf, R. C., Modern Control Systems, Addison-Wesley Publishing Company, 1992. 


66 


INITIAL DISTRIBUTION LIST 


Defense Technical Information Center 
Cameron Station 
Alexandria VA 22304-6145 


Library, Code 052 
Naval Postgraduate School 
Monterey CA 93943-5002 


Chairman, Code ME 

Department of Mechanical Enginnering 
Naval Postgraduate School 

Monterey CA 93943-5000 


Professor Fotis A. Papoulias, Code ME/PA 
Department of Mechanical Engineering 
Naval Postgraduate School 

Monterey CA 93943-5000 


Deniz Kuvvetleri Komutanligi 
Personel Daire Baskanligi 
Bakanliklar, Ankara/TURKEY 


Golcuk Tersanesi Komutanligi 
Golcuk, Kocaeli/TURKEY 


Deniz Harp Okulu Komutanligi 
81704 Tuzla, Istanbul/TURKEY 


Taskizak Tersanesi Komutanligi 
Kasimpasa, Istanbul/TURKEY 


Zeki O. Oral 


6. Cadde 35/12 
Bahcelievler, Ankara/TURKEY 


67 


No. Copies 
2 











DUDLEY KNOX LIBRARY 


NAVAL PGST GRADUATE SCHOOL 
MONTEREY CA 93943-5101 





JUDLEY KNOX LIBRARY 


’ i i } is | ae 
1 4 { phyrdy i wie { at ‘ted 2 La ae a - 
v; od at = 4 ni ite ( up E - , 1 ) a ét ie : lt tA! : 
rand itip assent ey Na: HYMN LIMA QUIN LINUEUT HU tt de Tt ! OR OHT: 
Oa Ee \ ‘ t | | | : 4, Nee ag i : 4 
hp ving So8 uN t S atts tale | l) \ | ‘ } , per. 
rey heast | seas ; ‘ j 
ea lf ith v pret : ae " | | | | | | | if } | | | | | ! * sank i tal Hae yaa + "ot hs a 
. bY a » | . £ : f 3 st a § 
i Lf d ; | | ] | | | | it ' thet aed vase 4 as a ty 8 ata Lv ehew 7 ey Aa atten’ a 
oaf 4 7 t | Vth | | | | | pg at J HPRaL 5 P Beye det er Maes Me ae F He 
a ti ae 2 Raaic AS yas | {| | | ia 1ibtil | | pea ats rh sone fapave pf bed ee aera 
; iebety! , rey i ; of bee, dee eee A Ye 
« H 2 : . ' 
4 ‘ ‘ 
oe H 1 ’ 
i, . 


2768 00308423 7 _ iam 


abe Rate Ms » ay it : - 
ry Artery > be Sk Uy ced bale ei pe t a pr peeeaded, beige int: : 
be aN ; * ad Py “4 ‘ . "hs ti : HP et Bes or 
- Loses al Hbboter ’ j : 1 
J the é wi : 
eget t ; 
yereisedee - 1 yt eeonet boa a's or a Ea 
2 ‘ oh mcirg’ eo ‘ of reataraca’ 
ait : g , 3 4 rr a eke Oberst ed big tebe bo paceate dt Fu c 
M . . , . ae eee 
x ° vad * aye Al 4 bewed ghes fe t 
.e br ‘ . t tgh eRag. ro wephPe ant 
- fin fas bop ao © Beene e fe, 
or yi Stee £6 . 
wate Ese cP ed gp gah, sae aus veg greed : 
> # ogeb g.0nhe fine 1 nee Pe 
Be iis a ats a6 
iP ygetet beam exe? J 
yeah oy f 4 rytenebe side 2 Bait 
ae er er “ : lis | er biey *y weep yi wld cneg ‘ 
bapaaek rae se BME late, pied PTET tng, atts ite tieneteeet eae ) 
. a ‘ 
’ be yee i iw " } eek ideal acnntta § : f Ph ‘ ry = ; LB 1 
“tb ean ceeds ae A ’ If tn . eyate & 
. ‘ ea 
a 


ob art 
ot { b , : : ta e te ib 
UF rte 1 Ot Paest oe Orb fy ' z ? : esl bs 
e Pak Migs bide be ; 4 Ms 16 ae adele 
a he a ; v : ‘ 5 baat rier 
yk eee h ee é SAE aes Be pel Sad as? A ye ra u set eee astaabetentis , 
tA id weet Ce ott eal Pires 3) ey ea AL 
‘ ir. Ta pee ages te deh 


. yal ae 
tr’ see Spt the 


A sot 

‘ fie! ¢ ; rh + rh yi 61 ets Lyi 

ink be bf asihiery as * 4 soe Pers ens 

a : H a Syregrt] rs 

3 ght tx 
mgd 

age toa! tpt aha Sy foks } ob Bh 
ye uv. : tt he A 

Sretal wl paca a hae 

“F A Lae Fra8 4. et4 

o . onvelerecgsdathepics 


. 
fe 
iat 
t's vf ‘t 


ol phar? 
feet aby pmeted ant. torts ane ae} VB 
ing PEAT eo : Wn rei ¢ acai Live 4 ae 2 dont’, ‘ 
D ta ligrestrl- pur bu babyh eed hd % IY 1” ‘ 
“He t 4 A peetangad & ba boned a at Goi) 
irre o8 : epadtet ates & t h rH ‘ a fee Wega Feed yt 
ons ie , bt PepRl yd yeh Pe és “¢ fedett oe 
4a bet be oak be oe I 4 : peyabsen epee cork apetad py 44 A 
Sr ba doed ' re . 4 he OMCs a ae ‘ 3} v dpestrerd pre 
Forpayeaee ff Sy ea hae — ty i it Y ibys Pda € iH jyiteceviaaey r ‘, fs i 
suf tas 4 gir erst dined deuse’, 8 OF : J Ud 
a ATE A if a ec deehetias ad i Garde tad : ryt abea, tea! a : DUES haiti sear Le 
‘ Sue y “ee 2 ony @) pofutag ot 
vat ‘“ oaks ' ap : yas op tingatet “ 4 LL A rete 
} 1 : ba gket eat h - E 
uae astute Gey peat pit 1% De Aa ena Pr 
geP Ket cc ’ ‘ Mears ae as : A! ate eee at 
ite : ie eee gh a ee 4 \ it gre : ue ate Hen i th ay vis tybee 
Ay ad ue} I ; bile's 2D H aa pt are Be ane arees 
ab BP hee B : : t ; cd Fee an cysnhel 408: = ‘ ofa 
L H vifeves) green teres : nI70 ‘ 


4 
i 
: 


anager 
yeh dyriseel ote 
roorarhh he 
un s 
std 
PLPe sees ebee 
agent ola ate, P 
Mioth =f : f 
C theta 
ace ‘ with 


ee are 


Fogel 
a 


as 
gtr oteF 


he gers 
ae hsb eb ess hagas! 
; r H H H on hity. an) 
oad deh pg Sil oe | fey 4 4 4 eae ne t nae fF «Be foots E 
nt ay : t j 4 qi ats ' aapi ew Ont 
: r peepee tts re bar shh ;' ae ge 
é rat r : : PCE mp bap ; OH Io gubea 
$ : the +e 


peraly ne ¢ reel 
+ A Pas nhese 2 ae 

sepde Fedeiey * 

hes Fie sf : : 

shes : : to 2 id § arte gelt a 

- ¢ 4 . 
ines’. Nees DeEPas ‘ " mae 

Sedipyed tha cbtes . 3° ; , ei teg sett 

oh gagee als dey d i 

va 


oom 
7 oy 


ifs s 
“Fe pre ot ‘ 
ateal oo 
qema (CB pe 0? 
eh rts @ ALY % 
is ye 
1 9p 08 ‘ 
dws’ af tye! 


Sat 
ae 
Mae 
cr * neem ee = 


=e wae 
a8 Que ler he Pe me 


‘ 
oe endl tt. 


en ee 
>= 
“on 


r Met 
cubes 39 Sa} Br wil a tekee ther 
rad gseta Ph Mews tty . ; :  aabse 
ab yrgdraa, ibevA 1, , ae re ea 
Paty higr ttt ye ; ' ergoicas gad fy atk te 
ane = 5 Ermase brvreds eds, ¢ 
+ e ts 
i Tphynt fs, Pee e ey? 
1 


Sle, 


Hers as 


herd 


ae 


7am *§ °*h 
3 Ther ade? 
8” ’ 


ad aihkye 
Lee 


ye Se 


aus 
=o ew 


ce 


Paw * 


Ste OS VET “At mm 
2g to we =e 


a 
=e See tn 
- 
Tes 
= 
a 


raf fhe 

pe ’ 
eet da tyme 

= ie a” mg da *8 

tsp oy hey dart 


= hm 5 
= 


eS 
os 


ery Parnar se? 
‘ rst 


°e) ATM 
SU) panels 
Cte Dif oP 
fay t a 
PB ars 
GEc* Op whe tates 
aegis 


3: 
ay 


€ 
4 


avint Sov 
ate eS: “ 
Pratatiri ae ¥ . . s 3 : y 
3 4 i : - Straus) 
Poel) 
+ f > 
See vas Ry z ey orl Rees eae 
i . t, tire ad # geet 
gee eee 


paccte Gasinee 


o tease oF 
wpete Ft 
263 


A 


mp4 A ea 


etcn SS ue 


yas 
t ie “i 304 Lael > 
4 BU latad add | 
ri Senr. ° 6 ™ 
LaRue ies Oe at 

Sah 


ene \ 

tet Ae 

dent OPH 
< 7 
a? 


ays? ! . ‘ f atit»s 
: : ; t gabe pedi re, ose st 
, Tikes g°S 9 tity cat a sate’ : pened a fy 
A a aire t "% Hae ? ‘ + | at wv ponga ft 
- = % i J 
i 3 


vareyie gtd ev 
o*es . ® ° 
i ha et Fert acre eer Meret or 7) 
ot? bogrtee § +. "9 sa tstef 5 
ye aga , f SiSsaree ot Page 
; CLT Sar “I, be 
stupor gtegt 


‘oar 
“Ue Tape qpuce 
tg ce “patete Patan 
ae, etete 86 2% ta %4 
en 4s Pie ocd | 
Pe Oe Ok a 


‘ . 
tee stm ta bh M8 
¢ 


awrnekoa ax 


“= 


‘ 


cis es 


os 
pees 


SPabsew 
t 


f , 
Satay tat sive HP 
ei2-07* es tee ‘ a Ad F A : ; ‘ : ; 
tory 4 pes SAgS y : f F ree iar? 1 fie, 
j . 5 “4 
¥, car ayes} ye ATA een 
: , ane ee) ee , 
Pere es Pepes 
rat . 
sage pyl oi Fea ataay 
by +8. ita he 


rah ton 
aw 


z 
ee 


wease Tabet y 

pa Pade ‘ aya 

Pa) hae ore the r] . a aS 

ee fe by a one | 4 

Begiyeaeh et “ane a { “pret 

gape tee ft cb rae Pani ben) 
si tt) ed 

To nmagh rasa 

get h 


on ns?e 2 
ney 3. ' to : Pe i 
ns . 
rH volnty ty Le ees » ' Ht af SSshee se i j ‘ Pa 
\ » , 4 = i : 5 af 7%» P 
A iY : ' ; 9th yaad (a ‘ $ a Me Z a" 9°¢ 
4 +" taped - ° ry 2 a) Aa ’ a 
. é F , re OY Oya tH DRS 
y i sinrfts, R t wad a 


ee cergedonet 
aie as ee 
+7 

. ‘ 


+> 


% 


+ oe tP  wale 
i24@ ee 
LS | 
P trge€ co, 


Onn = 


Ms 
fe 
33° 7k 
ete tgtee 
H #y 


3 
Opre, hetee Fes Riss 
Se oe 


ta i 


at 
ites 
“Wy 





